%% term_Bold_rot_B_V
function matrix_A = term_Bold_rot_B_V(mesh, fespace, Bold)
Nb1 = fespace.Nb;
Nlb1 = fespace.Nlb;
rows = zeros(mesh.Ne*Nlb1*Nlb1,1);
cols = rows;
val1 = rows;
val2 = rows;
val3 = rows;
val4 = rows;
flag = 1;
for n = 1:mesh.Ne
    V = mesh.node(mesh.elem(n,:),:);
    [int_points_loc, int_weights_loc] = integral_loc(V, "2d");
    phi = basis_loc(V, int_points_loc(:,1), int_points_loc(:,2), fespace.basis_type, "x");
    phix = basis_loc(V, int_points_loc(:,1), int_points_loc(:,2), fespace.basis_type, "dx");
    phiy = basis_loc(V, int_points_loc(:,1), int_points_loc(:,2), fespace.basis_type, "dy");
    B1 = phi * Bold(fespace.elem(n,:));
    B2 = phi * Bold(fespace.elem(n,:)+Nb1);
    for i = 1:Nlb1
        for j = 1:Nlb1
            f1 = -B2.*phiy(:,j).*phi(:,i);
            f2 = B2.*phix(:,j).*phi(:,i);
            f3 = B1.*phiy(:,j).*phi(:,i);
            f4 = -B1.*phix(:,j).*phi(:,i);
            r1 = int_weights_loc * f1;
            r2 = int_weights_loc * f2;
            r3 = int_weights_loc * f3;
            r4 = int_weights_loc * f4;
            rows(flag) = fespace.elem(n,i);
            cols(flag) = fespace.elem(n,j);
            val1(flag) = r1;
            val2(flag) = r2;
            val3(flag) = r3;
            val4(flag) = r4;
            flag = flag + 1;
        end
    end
end
A1 = sparse(rows, cols, val1, Nb1, Nb1);
A2 = sparse(rows, cols, val2, Nb1, Nb1);
A3 = sparse(rows, cols, val3, Nb1, Nb1);
A4 = sparse(rows, cols, val4, Nb1, Nb1);
matrix_A = [A1, A2; A3, A4];
end